Ciência de Dados para Negócios: Big Data for Finance Project


Resumo

teste de futuro para as ações

Intro

escrever

Code
library(tidyverse)
library(tidymodels)
library(modeltime)
library(timetk)
library(purrr)
library(tidyquant)
library(tsibble)
library(prophet)
library(feasts)
library(fable)
library(fabletools)
library(lubridate)
library(tictoc)

Carregamos os dados:

Code
tickers <- c(
         "BRFS3.SA",
  "JBSS3.SA",
  "BEEF3.SA",
  "MRFG3.SA",
  "TSN",
  "HRL",
  "GIS"
)

Então baixo os dados via Yahoo!Finance:

Code
portfolioPrices <- NULL
  for ( Ticker in tickers )
    portfolioPrices <- cbind(
      portfolioPrices, 
      quantmod::getSymbols.yahoo(
        Ticker,
        from = "2019-01-01",
        auto.assign = FALSE
      )[,4]
    )

portfolioPrices <- portfolioPrices[apply(portfolioPrices, 1, function(x) all(!is.na(x))),]

colnames(portfolioPrices) <- c(
  "BRFS3.SA",
  "JBSS3.SA",
  "BEEF3.SA",
  "MRFG3.SA",
  "TSN",
  "HRL",
  "GIS"
)

# Visualizar com DT
#DT::datatable(tail(portfolioPrices), options = list(pageLength = 10, scrollX = TRUE)) 

Visualizando os dados dos nossos últimos retornos dos preços, temos:

Code
log_returns <- log(portfolioPrices) - log(lag(portfolioPrices))
log_returns <- na.omit(log_returns)
log_returns <- log_returns |> 
  timetk::tk_tbl(preserve_index = TRUE, rename_index = "date")

tail(log_returns)
Code
ln_returns <- log_returns

ln_returns |> as.data.frame() |>
  dplyr::mutate(
    time = seq_along( TSN )
  ) |> select(-date) |>
  tidyr::pivot_longer(
    !time,
    names_to = "Variables",
    values_to = "Value"  
      ) |>
  dplyr::group_by(Variables) |>
  timetk::plot_time_series(
    time,
    Value,
    .interactive = F, # Change for TRUE for better visualization
    .facet_ncol = 2,
    .smooth = FALSE
  ) +
  ggplot2::theme(
    strip.background = ggplot2::element_rect(fill = "white", colour = "white")
  )

Modelagem com fpp3 e validação cruzada temporal

Precisaremos fazer um forecasting de curto prazo com nossos dados históricos de retornos pra formularmos nossas recomendações posteriores de compra, venda e espera:

  • Vamos começar com uma série por vez \(\Rightarrow\) TSN
Code
# Primeiro converto pra tsibble

lnretTSN <- log_returns |> 
  select(date, TSN) |> 
  as_tsibble(index = date)

glimpse(lnretTSN)
Rows: 1,519
Columns: 2
$ date <date> 2019-01-03, 2019-01-04, 2019-01-07, 2019-01-08, 2019-01-09, 2019…
$ TSN  <dbl> 0.0211432803, 0.0122208208, 0.0160060857, 0.0264099860, -0.017528…
Code
treino <- lnretTSN |>
  filter_index(~"2025-01-01")

War models

Code
tic()

Modelos <- treino |>
  model(
    AjusteExp = ETS(TSN ~ error("A") + trend("N") + season("N")), # Ajuste Exponencial com auto
    
    AjExp_aditivo = ETS(TSN ~ error("A") + trend("A") + season("A")), # Ajuste Exponencial Aditivo
    
    AjExp_multiplicativo = ETS(TSN ~ error("M") + trend("A") + season("M")), # Ajuste Exponencial Multiplicativo
    
    Croston = CROSTON(TSN), # Modelo Croston
    
    HoltWinters = ETS(TSN ~ error("M") + trend("Ad") + season("M")), # Holt Winters
    
    Holt = ETS(TSN ~ error("A") + trend("A") + season("N")), # Holt
    
    HoltAmort = ETS(TSN ~ error("A") + trend("Ad", phi = 0.9) + season("N")), # Holt Amortecida
    
    Regr_Comp = TSLM(TSN ~ trend() + season()), # Regressao com tendencia e sazonalidade auto
    
    Regr_Harmonica = TSLM(TSN ~ trend() + fourier(K = 2)), # Regressao harmonica
    
    Regr_Quebras = TSLM(TSN ~ trend(knots = c(2018, 2019, 2020))), # Regressao com quebras estruturais
    
    Snaive = SNAIVE(TSN), # SNAIVE
    
    Naive = NAIVE(TSN), #NAIVE
    
    Media_Movel = ARIMA(TSN ~ pdq(0,0,1)), # Media Movel Simples
    
    autoARIMA = ARIMA(TSN, stepwise = FALSE, approx = FALSE), # Auto ARIMA
    
    autoARIMA_saz = ARIMA(TSN, stepwise = FALSE, approx = FALSE, seasonal = TRUE), # AutoARIMA Sazonal
    
    #    Regr_erros_ARIMA = auto.arima(TSN, xreg = fourier(K = 3), seasonal = FALSE), # Regressao com erros ARIMA
    
    ARIMA_saz_012011 = ARIMA(TSN ~ pdq(0,1,2) + PDQ(0,1,1)), # ARIMA Sazonal ordem 012011
    
    ARIMA_saz_210011 = ARIMA(TSN ~ pdq(2,1,0) + PDQ(0,1,1)), # ARIMA Sazonal ordem 210011
    
    ARIMA_saz_0301012 = ARIMA(TSN ~ 0 + pdq(3,0,1) + PDQ(0,1,2)), # ARIMA sazonal
    
    ARIMA_quad = ARIMA(TSN ~ I(trend()^2)), # ARIMA com tendencia temporal quadratica
    
    ARIMA_determ = ARIMA(TSN ~ 1 + trend() + pdq(d = 0)), # ARIMA com tendencia deterministica
    
    ARIMA_estocastico = ARIMA(TSN ~ pdq(d = 1)), # ARIMA com tendência estocastica
    
    Regr_Harm_dinamica = ARIMA(TSN ~ fourier(K=2) + PDQ(0,0,0)), # Regressao Harmonica Dinamica
    
    Regr_Harm_Din_MultSaz = ARIMA(TSN ~ PDQ(0, 0, 0) + pdq(d = 0) + fourier(period = 7*30, K = 10) + fourier(period = 7*30, K = 5)), 
    
    Regr_Harm_Din_Saz = ARIMA(TSN ~ PDQ(0, 0, 0) + pdq(d = 0) + fourier(period = "month", K = 10) +
                                fourier(period = "year", K = 2) ), # Rgr Harm Mult Saz Complexa
    
#    Auto_Prophet = prophet(TSN), # Auto prophet
    
#    Prophet_mult = prophet(TSN ~ season(period = "month", order = 2, type = "multiplicative")),
    
#    Prophet_aditivo = prophet(TSN ~ season(period = "month", order = 2, type = "additive")),
    
#    Prophet_geom = prophet(TSN ~ growth("geometric") + season(period = "month", order = 2, type = "multiplicative")),
    
#    Prophet_memo = prophet(TSN ~ growth("geometric") + season(period = "month", order = 5) +
#                             season(period = "year", order = 2, type = "multiplicative")),
    
    Modelo_VAR = VAR(TSN, ic = "bic"), # Vetor Autoregressivo 
    
    Random_Walk = RW(TSN ~ drift()), # Random Walk com drift
    
    Rede_Neural_AR = NNETAR(TSN, bootstrap =  TRUE)#, # Rede Neural com auto AR e bootstraping nos erros
    
    #    x11 = X_13ARIMA_SEATS(TSN ~ x11()) # X11 ARIMA Seats
    
  ) |>
  
  forecast(h = "24 months") # Horizonte de projecao para os proximos 30 dias apos corte no treino

toc()  
1.02 sec elapsed

Selecionamos o melhor modelo (1 fold de validação cruzada somente):

Code
Modelos |>
  accuracy(lnretTSN) |>
  arrange(RMSE) # Seleção da acuracia pelo menor RMSE para o conjunto de modelos

Gero um cenário com o modelo:

Code
fit <- lnretTSN |>
  model(
    Regr_Quebras = TSLM(TSN ~ trend(knots = c(2018, 2019, 2020))), # Regressao com quebras estruturais
  )

sim <- fit |> generate(h = 30, times = 5, bootstrap = TRUE)

Plotamos os forecasts com esse modelo pra três cenários distintos no futuro:

Code
lnretTSN |>
  filter_index("2025-01-01"~.) |>
  ggplot(aes(x = date)) +
  geom_line(aes(y = TSN)) +
  geom_line(aes(y = .sim, colour = as.factor(.rep)),
    data = sim) +
  labs(title="Valores projetados de retornos de preços de contratos futuros da TSN", y="$US" ) +
  guides(colour = "none")

Code
# Primeiro converto pra tsibble

lnretGIS <- log_returns |> 
  select(date, GIS) |> 
  as_tsibble(index = date)

glimpse(lnretGIS)
Rows: 1,519
Columns: 2
$ date <date> 2019-01-03, 2019-01-04, 2019-01-07, 2019-01-08, 2019-01-09, 2019…
$ GIS  <dbl> 0.0154921374, 0.0200387411, 0.0164387227, 0.0149567729, -0.016440…
Code
treino <- lnretGIS |>
  filter_index(~"2025-01-01")
Code
tic()

Modelos <- treino |>
  model(
    AjusteExp = ETS(GIS ~ error("A") + trend("N") + season("N")), # Ajuste Exponencial com auto
    
    AjExp_aditivo = ETS(GIS ~ error("A") + trend("A") + season("A")), # Ajuste Exponencial Aditivo
    
    AjExp_multiplicativo = ETS(GIS ~ error("M") + trend("A") + season("M")), # Ajuste Exponencial Multiplicativo
    
    Croston = CROSTON(GIS), # Modelo Croston
    
    HoltWinters = ETS(GIS ~ error("M") + trend("Ad") + season("M")), # Holt Winters
    
    Holt = ETS(GIS ~ error("A") + trend("A") + season("N")), # Holt
    
    HoltAmort = ETS(GIS ~ error("A") + trend("Ad", phi = 0.9) + season("N")), # Holt Amortecida
    
    Regr_Comp = TSLM(GIS ~ trend() + season()), # Regressao com tendencia e sazonalidade auto
    
    Regr_Harmonica = TSLM(GIS ~ trend() + fourier(K = 2)), # Regressao harmonica
    
    Regr_Quebras = TSLM(GIS ~ trend(knots = c(2018, 2019, 2020))), # Regressao com quebras estruturais
    
    Snaive = SNAIVE(GIS), # SNAIVE
    
    Naive = NAIVE(GIS), #NAIVE
    
    Media_Movel = ARIMA(GIS ~ pdq(0,0,1)), # Media Movel Simples
    
    autoARIMA = ARIMA(GIS, stepwise = FALSE, approx = FALSE), # Auto ARIMA
    
    autoARIMA_saz = ARIMA(GIS, stepwise = FALSE, approx = FALSE, seasonal = TRUE), # AutoARIMA Sazonal
    
    #    Regr_erros_ARIMA = auto.arima(TSN, xreg = fourier(K = 3), seasonal = FALSE), # Regressao com erros ARIMA
    
    ARIMA_saz_012011 = ARIMA(GIS ~ pdq(0,1,2) + PDQ(0,1,1)), # ARIMA Sazonal ordem 012011
    
    ARIMA_saz_210011 = ARIMA(GIS ~ pdq(2,1,0) + PDQ(0,1,1)), # ARIMA Sazonal ordem 210011
    
    ARIMA_saz_0301012 = ARIMA(GIS ~ 0 + pdq(3,0,1) + PDQ(0,1,2)), # ARIMA sazonal
    
    ARIMA_quad = ARIMA(GIS ~ I(trend()^2)), # ARIMA com tendencia temporal quadratica
    
    ARIMA_determ = ARIMA(GIS ~ 1 + trend() + pdq(d = 0)), # ARIMA com tendencia deterministica
    
    ARIMA_estocastico = ARIMA(GIS ~ pdq(d = 1)), # ARIMA com tendência estocastica
    
    Regr_Harm_dinamica = ARIMA(GIS ~ fourier(K=2) + PDQ(0,0,0)), # Regressao Harmonica Dinamica
    
    Regr_Harm_Din_MultSaz = ARIMA(GIS ~ PDQ(0, 0, 0) + pdq(d = 0) + fourier(period = 7*30, K = 10) + fourier(period = 7*30, K = 5)), 
    
    Regr_Harm_Din_Saz = ARIMA(GIS ~ PDQ(0, 0, 0) + pdq(d = 0) + fourier(period = "month", K = 10) +
                                fourier(period = "year", K = 2) ), # Rgr Harm Mult Saz Complexa
    
#    Auto_Prophet = prophet(TSN), # Auto prophet
    
#    Prophet_mult = prophet(TSN ~ season(period = "month", order = 2, type = "multiplicative")),
    
#    Prophet_aditivo = prophet(TSN ~ season(period = "month", order = 2, type = "additive")),
    
#    Prophet_geom = prophet(TSN ~ growth("geometric") + season(period = "month", order = 2, type = "multiplicative")),
    
#    Prophet_memo = prophet(TSN ~ growth("geometric") + season(period = "month", order = 5) +
#                             season(period = "year", order = 2, type = "multiplicative")),
    
    Modelo_VAR = VAR(GIS, ic = "bic"), # Vetor Autoregressivo 
    
    Random_Walk = RW(GIS ~ drift()), # Random Walk com drift
    
    Rede_Neural_AR = NNETAR(GIS, bootstrap =  TRUE)#, # Rede Neural com auto AR e bootstraping nos erros
    
    #    x11 = X_13ARIMA_SEATS(TSN ~ x11()) # X11 ARIMA Seats
    
  ) |>
  
  forecast(h = "24 months") # Horizonte de projecao para os proximos 30 dias apos corte no treino

toc()  
0.99 sec elapsed
Code
fit <- lnretGIS |>
  model(
    Regr_Quebras = TSLM(GIS ~ trend(knots = c(2018, 2019, 2020))), # Regressao com quebras estruturais
  )

sim <- fit |> generate(h = 30, times = 5, bootstrap = TRUE)
Code
lnretTSN |>
  filter_index("2025-01-01"~.) |>
  ggplot(aes(x = date)) +
  geom_line(aes(y = TSN)) +
  geom_line(aes(y = .sim, colour = as.factor(.rep)),
    data = sim) +
  labs(title="Valores projetados de retornos de preços de contratos futuros da SALESFORCE", y="$US" ) +
  guides(colour = "none")

 

 


References


Markowitz, H. (1952). Portfolio Selection. The Journal of Finance, 7(1), 77–91.
Link

Sharpe, W. F. (1966). Mutual Fund Performance. The Journal of Business, 39(1), 119–138.
Link

Elton, E. J., Gruber, M. J., Brown, S. J., & Goetzmann, W. N. (2007). Modern Portfolio Theory and Investment Analysis (9th ed.). Wiley.

Hilpisch, Y. (2018). Python for Finance: Mastering Data-Driven Finance. O’Reilly Media.